#include <assert.h>
#include <math.h>
#include <stdio.h>

#include "../src/multidata.c"

static void
fill (float*           x,
      double         (*func) (void),
      unsigned int      n)
{
  for (unsigned int i = 0; i < n; ++i)
    {
      x[i] = func ();
    }
}

static void
fill_i (int*            x,
        int           (*func) (void),
        unsigned int    n)
{
  for (unsigned int i = 0; i < n; ++i)
    {
      x[i] = func ();
    }
}

static void
fill_r (Rating*         x,
        unsigned int    n_workers,
        unsigned int    n_labels,
        unsigned int    n_images,
        unsigned int    n)
{
  for (unsigned int i = 0; i < n; ++i)
    {
      x[i].worker = (lrand48 () >> 16) % n_workers;
      x[i].image = (lrand48 () >> 16) % n_images;
      x[i].rating = (lrand48 () >> 16) % n_labels;
    }
}

static double
zero (void)
{
  return 0.0;
}

static double
one (void)
{
  return 1.0;
}

static int
minus_one (void)
{
  return -1.0;
}

static int
unsigned_int_cmp (int a,
                  int b)
{
  if (a < b) { return -1; }
  else if (a > b) { return 1; }
  else { return 0; }
}

static int
float_cmp (float a,
            float b)
{
  if (a < b) { return -1; }
  else if (a > b) { return 1; }
  else { return 0; }
}

static int
fuzzy_cmp (const float* a,
           const float* b,
           unsigned int  n)
{
  const float epsilon = 1e-5;

  for (unsigned int i = 0; i < n; ++i)
    {
      int rv = fabs (a[i] - b[i]) < epsilon ? 0 : float_cmp (a[i], b[i]);

      if (rv) { return rv; }
    }

  return 0;
}

static int
data_cmp (const MultiData* a,
          const MultiData* b)
{
  return  unsigned_int_cmp (a->raw_labels, b->raw_labels) ||
          unsigned_int_cmp (a->max_depth, b->max_depth) ||
          unsigned_int_cmp (a->base.n_workers, b->base.n_workers) ||
          unsigned_int_cmp (a->base.rank, b->base.rank) ||
          unsigned_int_cmp (a->base.n_labels, b->base.n_labels) ||
          unsigned_int_cmp (a->base.n_images, b->base.n_images) ||
          unsigned_int_cmp (a->base.n_ratings, b->base.n_ratings) ||
          float_cmp (a->base.priorloga.sq_stddev, b->base.priorloga.sq_stddev) ||
          float_cmp (a->base.priorgamma.mean, b->base.priorgamma.mean) ||
          float_cmp (a->base.priorgamma.sq_stddev, b->base.priorgamma.sq_stddev) ||
          float_cmp (a->base.priorlogbeta.mean, b->base.priorlogbeta.mean) ||
          float_cmp (a->base.priorlogbeta.sq_stddev, b->base.priorlogbeta.sq_stddev) ||
          fuzzy_cmp (&a->base.priorc.mean, &b->base.priorc.mean, 1) ||
          float_cmp (a->base.priorc.sq_stddev, b->base.priorc.sq_stddev) ||
          float_cmp (a->base.priord.mean, b->base.priord.mean) ||
          float_cmp (a->base.priord.sq_stddev, b->base.priord.sq_stddev) ||
    /*
           memcmp (a->base.loga,
                   b->base.loga,
                   a->base.n_workers * a->base.rank * sizeof (*a->base.loga)) ||
           memcmp (a->base.gamma, b->base.gamma, a->base.rank * sizeof (*a->base.gamma)) || 
    */
         memcmp (a->base.logbeta, b->base.logbeta, a->base.n_images * sizeof (*a->base.logbeta)) || 
    /*
         memcmp (a->base.c, b->base.c, a->base.n_labels * a->base.rank * sizeof (*a->base.c)) || 
         memcmp (a->base.d, b->base.d, a->base.n_labels * a->base.rank * sizeof (*a->base.d)) || 
    */
         unsigned_int_cmp (a->base.n_clamp, b->base.n_clamp) ||
         memcmp (a->base.clamp, b->base.clamp, a->base.n_images * sizeof (*a->base.clamp)) || 
         fuzzy_cmp (a->base.priorz, b->base.priorz, a->base.n_labels) ||
         memcmp (a->base.ratings, b->base.ratings, a->base.n_ratings * sizeof (*a->base.ratings)) ||
         memcmp (a->base.image_ratings,
                 b->base.image_ratings,
                 a->base.n_images * sizeof (*a->base.image_ratings)) ||
         memcmp (a->base.worker_ratings,
                 b->base.worker_ratings,
                 a->base.n_workers * sizeof (*a->base.worker_ratings));
}

static void
test_data_inverse_once (void)
{
  /* read_clamp (read_data (print_data (data), print_clamp (data))) == data */

  MultiData data;
  MultiData data_dup;

  data.raw_labels = 1 + (lrand48 () >> 16) % 9;
  data.priorlabel = malloc (data.raw_labels * sizeof (*data.priorlabel));
  for (unsigned int i = 0; i < data.raw_labels; ++i)
    {
      data.priorlabel[i] = drand48 ();
    }

  do
    {
      data.max_depth = 1 + (lrand48 () >> 16) % 9;
    }
  while (data.max_depth > data.raw_labels);
  data.base.n_labels = 0;
  for (unsigned int m = 0; m <= data.max_depth; ++m)
    {
      data.base.n_labels += binomial (data.raw_labels, m);
    }

  data.base.n_workers = 1 + (lrand48 () >> 16) % 9;
  data.base.rank = 1 + (lrand48 () >> 16) % 9;
  data.base.n_images = 1 + (lrand48 () >> 16) % 9;
  data.base.n_ratings = 1 + (lrand48 () >> 16) % 9;
  data.base.priorloga = parameterized_normal_distribution (1);
  data.base.priorgamma = normal_distribution (0, 1);
  data.base.priorlogbeta = normal_distribution (1, 1);
  data.base.priorc = normal_distribution (1.0 / data.base.rank, 1);
  data.base.priord = normal_distribution (1.0, 1);

  data.base.loga = malloc (data.base.n_workers * data.base.rank * sizeof (*data.base.loga));
  fill (data.base.loga, one, data.base.n_workers * data.base.rank);
  data.base.gamma = malloc (data.base.rank * sizeof (*data.base.gamma));
  fill (data.base.gamma, zero, data.base.rank);
  data.base.logbeta = malloc (data.base.n_images * sizeof (*data.base.logbeta));
  fill (data.base.logbeta, one, data.base.n_images);
  data.base.c = malloc (data.base.n_labels * data.base.rank * sizeof (*data.base.c));
  fill (data.base.c, one, data.base.n_labels * data.base.rank);
  data.base.d = malloc (data.base.n_labels * data.base.rank * sizeof (*data.base.d));
  fill (data.base.d, one, data.base.n_labels * data.base.rank);
  data.base.pz = malloc (data.base.n_images * data.base.n_labels * sizeof (*data.base.pz));
  fill (data.base.pz, drand48, data.base.n_images * data.base.n_labels);
  data.base.priorz = malloc (data.base.n_labels * sizeof (*data.base.priorz));
  for (unsigned int i = 0; i < data.base.n_labels; ++i)
    {
      unsigned int labels[data.raw_labels];
      unsigned int depth = index_to_set (i, labels, data.raw_labels);
      float logpz = 0.0;

      for (unsigned int j = 0; j < depth; ++j)
        {
          logpz += log (data.priorlabel[labels[j]]);
        }

      data.base.priorz[i] = logpz;
    }

  safe_exponentiate_normalize (data.base.priorz, data.base.n_labels);

  data.base.n_clamp = 0;
  data.base.clamp = malloc (data.base.n_images * sizeof (*data.base.priorz));
  fill_i (data.base.clamp, minus_one, data.base.n_images);
  for (unsigned int i = 0; i < data.base.n_images; ++i)
    {
      if (drand48 () < 0.5) 
        { 
          data.base.clamp[i] = (int) (drand48 () * data.base.n_labels); 
          ++data.base.n_clamp;
        }
    }

  data.base.ratings = malloc (data.base.n_ratings * sizeof (*data.base.ratings));
  fill_r (data.base.ratings,
          data.base.n_workers,
          data.base.n_labels,
          data.base.n_images,
          data.base.n_ratings);

  data.base.image_ratings = malloc (data.base.n_images * sizeof (*data.base.image_ratings));
  memset (data.base.image_ratings, 0, data.base.n_images * sizeof (*data.base.image_ratings));
  for (unsigned int i = 0; i < data.base.n_ratings; ++i)
    {
      ++data.base.image_ratings[data.base.ratings[i].image];
    }

  data.base.worker_ratings = malloc (data.base.n_workers * sizeof (*data.base.worker_ratings));
  memset (data.base.worker_ratings, 0, data.base.n_workers * sizeof (*data.base.worker_ratings));
  for (unsigned int i = 0; i < data.base.n_ratings; ++i)
    {
      ++data.base.worker_ratings[data.base.ratings[i].worker];
    }

  FILE* fp = tmpfile ();

  print_multi_data (fp, &data);
  print_multi_clamp (fp, &data);
  fflush (fp);
  rewind (fp);
  read_multi_data (fp, &data_dup, data.base.rank);
  read_multi_clamp (fp, &data_dup);

  assert (data_cmp (&data, &data_dup) == 0);

  fclose (fp);
  multi_data_destruct (&data_dup);
  multi_data_destruct (&data);
}

static void
test_data_inverse (void)
{
  for (unsigned int i = 0; i < 100; ++i)
    {
      test_data_inverse_once ();
    }
}

int
main (void)
{
  srand48 (69);
  test_data_inverse ();

  return 0;
}
